APS/123-QED 



Ensemble v-representable ab-initio density-functional calculation of energy and spin 
in atoms: a test of exchange-correlation approximations. 

Eli Kraisler 

Raymond and Beverley Sadder Faculty of Exact Sciences, 
School of Physics and Astronomy, Tel Aviv University, Tel Aviv 69978, Israel and 
Physics Department, NRCN, P.O.Box 9001, Beer Sheva 84190, Israel 

Guy Makov 

Department of Materials Engineering, Ben-Gurion University of the Negev, P.O.Box 653, Beer-Sheva 84105, Israel 

Itzhak Kelson 
Raymond and Beverley Sackler Faculty of Exact Sciences, 
School of Physics and Astronomy, Tel Aviv University, Tel Aviv 69978, Israel 
(Dated: November 23, 2010) 

The total energies and the spin states for atoms and their first ions with Z = 1 — 86 are calculated 
within the the local spin-density approximation (LSDA) and the generalized-gradient approxima- 
tion (GGA) to the exchange-correlation (xc) energy in density-functional theory. Atoms and ions 
for which the ground-state density is not pure-state w-representable, are treated as ensemble v- 
representable with fractional occupations of the Kohn-Sham system. A newly developed algorithm 
which searches over ensemble v-representable densities [E. Kraisler et al, Phys. Rev. A 80, 032115 
(2009)] is employed in calculations. It is found that for many atoms the ionization energies obtained 
with the GGA are only modestly improved with respect to experimental data, as compared to the 
LSDA. However, even in those groups of atoms where the improvement is systematic, there remains 
a non-negligible difference with respect to the experiment. The ab-initio electronic configuration in 
the Kohn-Sham reference system does not always equal the configuration obtained from the spec- 
troscopic term within the independent-electron approximation. It was shown that use of the latter 
configuration can prevent the energy-minimization process from converging to the global minimum, 
e.g. in lanthanides. The spin values calculated ab-initio fit the experiment for most atoms and are 
almost unaffected by the choice of the xc-functional. Among the systems with incorrectly obtained 
spin there exist some cases (e.g. V, Pt) for which the result is found to be stable with respect to small 
variations in the xc-approximation. These findings suggest a necessity for a significant modification 
of the exchange-correlation functional, probably of a non-local nature, to accurately describe such 
systems. 



PACS numbers: 31.15.A-, 31.15.ac, 31.15.E- 

I. INTRODUCTION. 

Density-functional theory (DFT) [U-IH is the leading 
theoretical framework for studying the electronic prop- 
erties of matter. Its approach is ab-initio, which means 
that in principle, the only necessary input for the theory 
is the external potential v(r) and the total number of elec- 
trons N; no experimental data is required. DFT is cur- 
rently applied to a variety of many-electron systems @,@] : 
along with traditional applications to atoms, molecules 
and solids, DFT is used for studying nano-objects, va- 
cancies and impurities, surfaces and even DNA. The ap- 
plications in all these fields are, obviously, interrelated. 

The many-body physics of Coulomb-interacting elec- 
tron systems is represented in DFT by the exchange- 
correlation (xc) energy functional E xc . It is non-local, 
spin-dependent, it has non-analytical properties, and its 
exact form is not known [1, d, Q. There exist, however, 
many approximations to this functional (e.g. (lp| - [l3 |). 
based on numerical results for the homogeneous electron 
gas [HI and asymptotic analytical derivations (see [l4| 
and references therein). The validity of alternative ap- 



proximations can only be determined by testing them on 
a wide range of systems, in comparison to the experimen- 
tal data. 

Several reasons motivate us to explore atomic systems 
with DFT. First, atoms are possibly the simplest physical 
systems DFT can describe, and they are relatively easy 
to solve numerically. Second, in atoms a wide spectrum 
of electron-electron interactions takes place, e.g. closed- 
shell, open-shell, magnetic, relativistic, strongly corre- 
lated systems and others. Such a diversity creates a chal- 
lenge to any approximation to the exchange-correlation 
energy. Third, there exists a large and highly accurate 
experimental database on atoms and ions |16| . Both the 
ionization energies and the spectroscopic terms can be 
found for all atoms in the periodic table. Therefore, a sys- 
tematic comparison of numerical findings to precise ex- 
perimental results is possible. For these reasons, atomic 
systems provide a " playground" to test physical approxi- 
mations and numerical algorithms within DFT. Last but 
not least, since atoms serve as building blocks of more 
complicate material systems, a deeper understanding of 
atoms and ions, and as a result - their better description, 
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can improve our ability to model real materials. 

Atomic systems have been extensively investigated 
with DFT in the last four decades. As early as 1966, 
Tong and Sham [l7j published density-functional calcu- 
lations for several alkali and rare gas atoms using a spin- 
unpolarized local density approximation. The study by 
Moruzzi, Janak and Williams [l8| (although mainly fo- 
cused on crystals) presented numerical results for a large 
group of atoms. The total energies of atoms were ob- 
tained within a spin-polarized version of DFT, using the 
LDA form by von Barth and Hedin |19]. In this work 
and in subsequent work by Janak |20| . the problem of 
fractional occupations in the Kohn-Sham reference sys- 
tem for the Fe and Co atoms was discussed explicitly. 
However, since the physical significance of unequal occu- 
pations of degenerate levels in the reference system was 
not fully understood at that time, the authors [3 H3] 
questioned their findings in this regard. Furthermore, 
these initial studies did not include any results for ions, 
nor the ionization energies of atoms. 

Subsequent to several partial studies [HI, [2lT[23 | , Ko- 
tochigova et al. (2f| analyzed all atoms and their first 
ions in the periodic table with local approximations to 
the xc-energy. These approximations included both spin- 
polarized and unpolarized functionals, within both rela- 
tivistic and non-relativistic regimes. The agreement of 
the LSDA ionization energies with experimental data was 
found to be reasonable. 

However, the electronic configurations of the atoms 
and ions employed in that study were not obtained ab- 
initio. For some systems the energy-minimization pro- 
cess was restricted in a manner that prevented it from 
converging to the global minimum and may produce er- 
roneous results, as will be discussed below. 

All results in the aforementioned works were obtained 
using the local (spin-) density approximation (L(S)DA). 
A systematic investigation of atomic systems with func- 
tionals constructed using the generalized gradient ap- 
proximation (GGA) has not yet been published. How- 
ever, Lee and Martin (2(| have calculated the total en- 
ergies for atoms with Z = 1 — 20, 31 — 36 obtained in 
an unpolarized calculation with PBE [3 and PW91 pj] 
GGA's and with LSDA(PZ) approximations [il) . Ion- 
ization energies for atoms with Z = 1 — 20 in a spin- 
polarized calculation were also presented. The authors 
found that PBE-GGA results for the total energies im- 
prove the LDA results considerably, while for the ioniza- 
tion energies PBE-GGA is only slightly better than the 
LSDA. 

Surprisingly, none of the studies mentioned above, re- 
gardless of the approximations used, analyzes how well 
do DFT calculations predict the spin of atoms and ions. 

Since the latest comprehensive works on density- 
functional theory of atoms, two major advances have 
taken place. First, questions of v-representability and 
differentiability of density-functionals in DFT have been 
clarified. With respect to atomic calculations, (a) the use 
of unequal fractional occupation of degenerate energy lev- 



els in the Kohn-Sham reference system has been found to 
be mathematically allowed and physically justified [9I.I27I- 
30]; (b) the use of excited states of the Kohn-Sham sys- 
tem in order to describe a ground state of the real sys- 
tem has been shown to be forbidden [301 ] . Second, a nu- 
merical algorithm designed to treat atomic systems with 
ensemble w-representable ground-state densities was pro- 
posed [30(, which allows for fractional occupation of de- 
generate KS levels. The algorithm was illustrated for the 
Fe atom - a well-known ensemble ii-representable (EVR) 
system [H, H3, HH . 

To conclude, the limitations of previous work and re- 
cent developments in density-functional theory suggest 
the need for a new comprehensive study for atoms and 
ions with DFT, which will ab-initio determine the KS 
electronic configuration, making it possible to obtain the 
spin, the total and ionization energies without relying on 
experimental findings, thus testing the approximations 
to the xc-energy employed. 

In the current contribution we present results of self- 
consistent ab-initio non-relativistic DFT calculations for 
atoms and their first ions with 1 Z 86, within the 
spherical approximation for the density. The objectives 
of this paper are: (i) to obtain ab-initio the spin, S, of 
the atoms and ions in the Kohn-Sham DFT scheme and 
compare it with experimental results; (ii) to explore the 
difference between the total and the ionization energies 
obtained while determining the occupation numbers ab- 
initio to those obtained using the empirically assigned 
occupation numbers [251 ]; (iii) to systematically test the 
effect of including the PBE-GGA on the energies and 
spins of atoms and ions, relative to the LSDA; (iv) to 
examine the quantitative difference between the ioniza- 
tion energies and the spins obtained while using various 
parametrizations for the correlation energy |10h12| . 

The rest of the paper is organized as follows. Sec. |TT] 
provides the theoretical background, in Sec. |HT] the nu- 
merical details concerning atomic calculations are dis- 
cussed and Sec. IIVI presents the numerical results for all 
atoms and ions with 1 ^ Z ^ 86. These include the 
total energy E, the first ionization energy /, the spin S 
and the Kohn-Sham electronic configuration. The results 
obtained are c omp ared to previous calculations from the 
literature [25l . |26| . as well as to experimental measure- 
ments [IH . Sec. IVl contains a discussion. Possible rea- 
sons for existing discrepancies are provided, and ways to 
reduce those in future research are proposed. 



II. THEORY. 

A non-relativistic many-electron system is de- 
scribed 043 by the Hamiltonian 
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where m e is the electron mass, N is the number of elec- 
trons in the system, v(f) is the external potential, which 
equals —Ze 2 /r for atoms and ions. 

According to the Hohenberg-Kohn (HK) theorems [lj, 
IH , the ground state density determines all the properties 
of a many-electron system. In the spin-dependent ver- 
sion of DFT, the total-energy functional E v [rif , nj.] can 
be represented as 



E v [rif, Hi] = / v n d 3 r + Fhk [%, n{\ , 



(2) 



where and n± are the partial densities of electrons with 
the spin a = t and I, n = iVf + rij. is the total density 
and Fhk is the Hohenberg-Kohn functional, which is in- 
dependent of the particular electron system, i.e. of v. An 
explicit expression for Fhk is not known, and it is usually 
approximated using the Kohn-Sham (KS) approach 
introducing an auxiliary system of non-interacting elec- 
trons subject to an effective potential v e jj, which is cho- 
sen so that the density of the KS system equals the den- 
sity of the interacting system. In spin-DFT, two coupled 
auxiliary systems are used to reproduce n-j- and n±. 

It is well-known that the densities in some atomic 
systems are ensemble- (rather than pure-state-) v- 
representable [H, [2(| El • Therefore, in the follow- 
ing all the expressions are generalized to include ensem- 
ble s, emp loying the results of previous studies 0, [27l — 
l30l . I3214341 ] . In particular, Lieb's extension to the HK 
functional [l, [27], HH is used, and it is denoted Fl below. 

Within the KS system, solving one-particle 
Schrodinger equations 



(3) 
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one obtains the eigenvalues {£ia} and the orbitals {ipia} 
of the KS system, some of which can be degenerate. 
The total energy of the c-KS system is given by 



E K s [riff] 



v e ff,anad 3 r + T L [na}, 



(4) 



where Tl is Lieb's kinetic energy functional 0, H3, HH 
for the KS system. 

The expression for the effective potential is found by 
applying Eulcr's differential relations to both the inter- 
acting and the KS systems: 
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The effective potential has the form (see e.g. [J]) 
v effA n t^ n l\( r l = v (r) + v H [n] + Vacant, rij, where 
VH[n] = e 2 J n(r')/\r — r'\d 3 r' is known as the Hartree 
potential and 
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is the exchange-correlation potential. The total energy of 
the interacting system can be expressed in terms of the 
KS system as 

E = T L [n t ] + T L [nj + J vnd 3 r + E H [n} + E xc [n t , , 

(8) 

where Eh and E xc are the Hartree and the exchange- 
correlation energies, respectively. 

It is evident, therefore, that in the KS approach the 
functionals Fl and Tl have to be differentiable with re- 
spect to the densities n a . The question of differentiability 
was discussed extensively in the literature [9|, I27H301 |32| — 
H- 

It was found that the partial densities must be non- 
interacting ensemble w-representable (NI-EVR), i.e. to 
have the form 



iff I > 



(9) 



where gi a are the occupation numbers, which obey 
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Ejff < £F<t 
&ia — £F<j 
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Here £i a is the energy of the ith KS level of the cr-system, 
Diu is the maximal number of electrons that can occupy 
the i-th level, Efo- and Xi a € [0,Di a ] are the energy and 
the occupation of the highest occupied level(s), which can 
be fractional. The kinetic energy then equals 



2l[«<t] 



2m, 
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It is emphasized [36] that the occupation numbers 
{gia-} are not free parameters, as was previously sug- 
gested 20]. All systems that satisfy this restriction are 
termed below proper. 

S, the z-projcction of the total spin (referred to in this 
article as spin) is a functional of the partial densities 

S=\J («t - H) d 3 r = I £> it - 9il ) = l(JV t - Nj. 

(12) 

Therefore, the spin of the KS system is identical to the 
spin of the interacting system. However, the occupation 
numbers {gia} are internal quantities of the KS systems 
determined by Eq. (flTJ|) . They are not necessarily equal to 
the occupation numbers that are available in the experi- 
mental literature for atoms [l6| ■ The latter (referred to in 
this work as empirical) are derived from the independent- 
electron model of the atom (see e.g. [jjj ch. X]) to satisfy 
the spectroscopic terms deduced from measured emission 
spectra, and cannot be explicitly related to their coun- 
terparts in the KS system. 

The spin S is determined in the KS scheme by find- 
ing the minimum total energy for different spin values 
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of the non-interacting reference system. These calcula- 
tions are performed employing an approximate exchange 
correlation functional which is a possible source of error 
in determining the energy and therefore creates uncer- 
tainty in the calculated total spin. Since the spin is a 
half-integral number, which can increase or decrease by 
unit value only, then the validity of the calculated result 
is determined by the magnitude of the energy difference 
between the system with spin S and the same system 
with spin S ± 1. The larger this energy difference is, the 
less likely it leads to an incorrect conclusion regarding the 
spin of the system. If the energy difference is small, then 
a small change in the exchange correlation energy could 
cause the order of the energies to change and therefore 
the predicted total spin of the system (see Sees. IIV El and 
HVFl) . 

With respect to differentiability, two results were re- 
cently obtained [3(J. First, it was shown that the use 
of excited states of the KS system to describe a ground 
state of the interacting system is not allowed. It can be 
claimed that a density that emerges in some KS systems 
from an excited state (i.e. that does not obey Eq. ([TU]) ) 
can be obtained also as a ground state of another KS sys- 
tem, so the differentiability is assured. Then, however, 
Euler's relations ([S]),© cannot be satisfied for both KS 
systems simultaneously, so the Kohn-Sham method can- 
not be applied. Second, for spin-polarized DFT, each of 
the partial densities rif, nj, has to obey Eqs. (|9j) . (TTOj) . but 
there can be gaps in the occupation of the whole system 
if the first vacant level of one of the subsystems lies below 
the last occupied level of the another one. We call such 
systems proper in a broad sense (b.s.). This situation 
arises because we are constrained from having partial oc- 
cupation of states with a different spin, by the physical 
requirement that the spin is half-integer. 

It is emphasized here that a b.s. -proper density obeys 
all the required restrictions related to a rigorous defini- 
tion and differentiability of energy functionals, thus it 
can serve as a legitimate candidate for a solution of a 
many-electron system. A b.s. solution can be converted 
to a proper solution by shifting one of the effective poten- 
tials, say v e jj ^, by a constant, employing the ambiguity 
in definition of potentials. This transformation will not 
affect the densities, but will change the relative position 
of the t and 4- KS energy levels. However, such a shift 
is analogous, in a known sense, to applying a uniform 
magnetic field on the system. We confine ourselves in 
the current work to zero magnetic field, which makes us 
distinguish b.s. proper solutions from fully proper ones. 

In the current work the physical quantities obtained in 
DFT calculations are the partial densities nf and and 
the total energy E. The ionization energy 

I = E + -E, (13) 

and the spin S (see Eq. (TT2]) ) are derived from those quan- 
tities. Here E is the total energy of an atom and E + 
is that of the corresponding ion. These are opposed to 
quantities like v e ff a , tpi a and g^, which belong to 



the KS system and do not receive here a direct physical 
interpretation. Note, however, that although gi a belong 
to the latter, the sums X^(ffit ^ 54 ) have an interpreta- 
tion, being equal to N and 2S, correspondingly. 



III. NUMERICAL METHODS. 

As in earlier studies of atomic systems [HI, [5^, |38[ , in 
the current work the density was approximated by its 
spherical average n a (r) = (47r) _1 J Q n n a (r)dfl. Since v 
is spherical as well, the effective potentials are spheri- 
cal, too. This approximation reduces a three-dimensional 
problem to a one-dimensional problem. With a spherical 
potential, the wave functions take the form i^nimoif) = 
Rnia{r)Yi m (9 , cf>) , where Yi m (8,(f>) are the spherical har- 
monics and Rnitj{r) are radial wave functions that obey 
the following differential equation: 

// 2 / ( m e, , NN 1(1 + 1)\ 

Rnla + -Rnla+[ 2 ^2 ( £nl <? ~ V eff A r )) ^ J R nlv = 0. 

(14) 

Since the energy levels do not depend on m, D n i a = 2Z+1. 

To determine the occupation numbers Qni„, the al- 
gorithm AEVRA presented in Ref. [30] was used. 
In this algorithm, we start with a reasonable guess for 
v e ff a (r) and occupy the energy levels of the reference 
systems properly, given the total spin value S. Therefore, 
we are assured that the density generated by the refer- 
ence system obeys the known mathematical restrictions 
(see Sec. HI)) . This allows us to calculate the effective 
potential from this density. If we employ this potential 
directly to calculate a new density in EVR cases, then as 
predicted in (39[, the process starts alternating between 
two electronic configurations, and no convergence takes 
place. 

Instead, as the first step, we reduce the linear mix- 
ing coefficient for the effective potential each time the 
electronic configuration changes. By repeating this pro- 
cedure, it is possible to find a v e jj, with energy levels 
that coincide to any required accuracy A e , thus ensuring 
we enter the domain of NI-EVR densities. However, this 
solution is not necessarily self-consistent, and additional 
iterations in the domain of NI-EVR densities are required 
to obtain the final result. 

The second step in the algorithm continues the search 
in the subspace of potentials that generate NI-EVR den- 
sities while relaxing the constraint on integer occupation. 
Once an effective potential with two degenerate levels is 
found, it can produce a range of proper NI-EVR densi- 
ties by fractionally occupying the degenerate levels. We 
choose the occupations in a way that preserves the de- 
generacy of the energy levels for the next iteration [30j . 
but without constraining their values. 

In principle, the described process has to be performed 
separately for all values of the total spin S, finally choos- 
ing the one with the lowest energy. In practice, however, 
very high values of S are not expected to minimize the 
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energy of the system, if they produce a bound state at 
all. Therefore, in the current work, for atoms and ions 
within the p-block only spins obtained by all possible oc- 
cupations of the highest and p^ levels were considered; 
for d-block systems - only spins obtained by occupying 
the highest dJ, dr, s^, and levels, and for the /-block 
- by occupying f', f+, dJ, dr, s' and levels. 

In all calculations, a high convergence of 2 /ihartree for 
the total energy was obtained. For cases with a degener- 
ate ground state of the reference system, the degenerate 
energy levels were kept within A e = 0.1 mhartree of each 
other, and thus fractional occupation numbers have an 
estimated precision of 10~ 3 . To assure the desired accu- 
racy in energy, Eq. (|14j) was solved on a logarithmic grid 
with 16,000 points, on the interval (e~~ a /Z,L) in Bohr 
radii, with a = 13 and L — 35 for the LSDA and 45 for 
the GGA. 



xc S E (hartree) configuration proper 

Co LSDA 3/2 -1380.194503 [Ar]3d2. 849 4si. 151 yes 

Co LSDA 1/2 -1380.181419 [Ar]3dl4s? bs 

Co GGA 3/2 -1382.509230 [Ar]3d|. 832 4sJ. 168 yes 

Co GGA 1/2 -1382.495184 [Ar]3ri|4si bs 

Co+ LSDA 2 -1379.854391 [Ar]3da4si bs 

Co+ LSDA 1 -1379.896444 [Ar]344sJ] yes 

Co+ LSDA -1379.865191 [Ar]3d|4s!] yes 

Co+ GGA 2 -1382.182683 [Ar]3dl4sJ bs 

Co+ GGA 1 -1382.218966 [Ar]3d|4s!5 yes 

Co+ GGA -1382.185990 [Ar]3d|4sg yes 



TABLE I: The total energy and the KS electronic configu- 
ration for the Co atom and its first ion obtained within the 
LSDA and the GGA (see definitions in text). The cases are 
identified as proper or proper in a broad sense (bs). 



A. Co - an Example of an Atomic DFT 
Calculation. 



IV. RESULTS. 

In the course of the current work the total and ion- 
ization energies, the spin and the KS electronic con- 
figurations for atoms and first ions with 1 ^ Z ^ 86 
were obtained within the local spin-density approxima- 
tion (LSDA) and the generalized gradient approximation 
(GGA) by Perdew, Burke and Ernzerhof [14j . Calcula- 
tions for all atoms and ions were performed using the 
parametrization by Vosko, Wilk and Nusair (VWN) [l(| 
for the correlation energy of a homogeneous electron gas. 
Atoms and ions with 1 ^ Z ^ 56 or 71 ^ Z ^ 86 
(excluding the lanthanides) were also solved with two 
other widely-used parametrizations - the one by Perdew 
and Zunger (PZ) and the one by Perdew and Wang 
(PW92) [12| ■ These results are summarized in a database 
of atomic DFT calculations 40] . 

For the first time, extensive results for atoms and ions 
within the PBE-GGA [l4| are presented. In addition, the 
spin of atoms and ions is obtained ab-initio due to the 
use of AEVRA algorithm in calculations. 

The results are presented as follows. In Sec. IIV Al thc 
calculations for the Co atom and its first ion are presented 
in detail to explain how the atomic database was 
constructed. Sec. IIVBI surveys results of the ionization 
energy for atoms from all over the periodic table. The 
total energies for light atoms (1 ^ Z ^ 29) are also 
included and compared to experiment. In Sec. IIV Cl the 
spin values obtained ab-initio are reported for various 
atoms and ions and Sec. IIV Dl shows what is the effect of 
restricting the KS occupation numbers to their empirical 
values. Sec. IIV El displays the improvement of PBE-GGA 
over the LSDA and Sec. IIV Fl compares the numerical 
results obtained with different parametrizations for the 
correlation energy of a homogeneous electron gas. 



The atom of Co has Z = 27 electrons. Experimen- 
tally [l6j], its first ionization energy I exp — 7.88101 eV, 
which equals 0.289622 hartree. An experimental value for 
the total energy of an atom can be obtained by summa- 
tion of all its ionization energies, until a full ionization. 
For Co, E exp — —1393.4 hartree. The measured spin 
value for Co is S = 3/2 with [Ar]3d|4s} as its empirical 
electronic configuration [48j. For the ion Co + , 5=1 
with a configuration [ArjSd^SQ. 

According to Sec. Mil the values of S — |, h were con- 
sidered for Co and for Co+ - S = 2,1,0. The results 
are given in table U for the LSDA and the PBE-GGA, in 
the VWN parametrization [10( . The calculated total en- 
ergy fits the experiment within a high accuracy of 0.9% 
for the LSDA and 0.7% for the GGA. The spin is repro- 
duced correctly for both the atom and the ion in both 
approximations. For the ion the empirical and the KS 
electronic configurations coincide. However, for the Co 
atom the KS configuration is different from the empiri- 
cal one, and includes fractional occupation numbers. Co 
is one of many NI-EVR cases where the assumption that 
the empirical electronic configuration can be used in DFT 
calculations (see e.g. (25[) is shown to be unjustified. 

Table U demonstrates that an emergence of a NI-EVR 
solution for Co is not an artifact related to the local- 
density approximation for the exchange-correlation func- 
tional as speculated in Refs. [H, Hfll) because it appears 
also for the generalized gradient approximation, with 
slightly different occupations of the degenerate levels. 

From table U the first ionization energy can be calcu- 
lated for Co: it equals 0.298059 hartree for the LSDA 
and 0.290264 hartree for the GGA. Thus, employment 
of a gradient-dependent xc-functional in Co significantly 
improves the result for the ionization energy, reducing 
the relative error from 3% to 0.2%. 
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B. Total and Ionization Energies Throughout the 
Periodic Table - a Survey. 

The available experimental data fl6l | on ionization en- 
ergies allows comparison of the total energy E ca i c ob- 
tained in calculations to that measured experimentally, 
E eX p, for atoms with 1 ^ Z 29. As mentioned in 
Sec. I IV A[ E exp is obtained by summation of all the ion- 
ization energies of a given atom. The relative difference 



\AE/E\ — \(E calc — E exp )/E e 



(15) 



is plotted in Fig. [TJ as a function of Z, for the LSDA 
and the PBE-GGA [14]. In addition, numerical results 
obtained with the Hartree-Fock (HF) method [4l[ are 
given for comparison. From this figure one learns that 
the LSDA provides an accuracy < 1% for all the ex- 
amined atoms, except the lightest ones - 1 ^ Z ^ 5, 
where discrepancies can reach 4.5%. This can be ex- 
plained by the fact that LSDA, based on results for a ho- 
mogeneous electron gas, cannot describe accurately sys- 
tems with few electrons. Including a better description 
of inhomogeneous many-electron systems, as done in the 
GGA, improves over the LSDA results systematically and 
considerably, especially for light atoms. The relative im- 
provement of the GGA diminishes with increase in Z. It 
is also clearly seen that the LSDA has a slightly worse 
accuracy compared to HF, and that the GGA performs 
better than HF and LSDA. 

For both LSDA and GGA, we note that the relative 
difference tends to grow for Z > 15. It is probable that 
this error in the total energy is due to relativistic con- 
tributions which become significant for the inner core 
electrons at these values of Z. This assumption is sup- 
ported by the following elementary argument. Since the 
energy of the innermost electron can be approximated by 
Ei = — Z 2 e 2 / (2ao) , where ao is Bohr's radius, the ratio 
£i/(m e c 2 ) increases from 0.003% for Z = 1 to 2% for 
Z = 25. It is also in agreement with the published re- 
sults [25| for the total energy of atoms calculated with 
relativistic LDA. Nonetheless, I is determined by the 
outermost electron, relativistic corrections become sig- 
nificant only at large values of Z. 

Figure [2] shows the ionization energy / as a function 
of the atomic number Z for the LSDA(VWN) and the 
PBE-GGA(VWN). Experimental results [l6[ are shown 
for comparison. Figure [3] plots the relative errors 



calc 



L exp | / *exp 



(16) 



of the LSDA and GGA results with respect to the ex- 
periment. Here I ca i c stands for the calculated ionization 
energy and I exp - for the experimental result. Table |TT] 
summarizes the average relative errors of the ionization 
energies for different groups of atoms. 

From the figures above we see that many trends in 
the experimental I(Z) are reproduced by the calculated 
results, e.g. the sharp decreases after an electron shell 
is filled (Z — 11,19,37,55) and the more shallow dips 
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FIG. 1: (Color online). The relative error in the total energies 
of different atoms with respect to the experiment. Rhombi 
used for LSDA, circles - for GGA calculations (this work); 
squares - for HF method (Ref. HH). 
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TABLE II: The average relative errors A re ; of the ionization 
energies obtained within the LSDA and the PBE-GGA with 
respect to the experiment. 



when the p-shell is half-full (Z = 8, 16, 34, 52), which are 
a straightforward consequence of a spin-polarized treat- 
ment. 

Table HT1 shows that ionization energies are reproduced 
within a few percent. For heavy atoms, the average accu- 
racy drops, however, and reaches 10%, due to relativistic 
effects. 

In lanthanides the results for the ionization energy are 
quite surprising. Despite the relatively simple form of 
the experimental curve, consisting of a straight line with 
a single small peak for Gd (Z = 64), the calculations yield 
a completely different shape, and the relative errors can 
reach 17% (see Sees. II V Dl and IVl for more details). 

To summarize, GGA improves over the LSDA for s-, p- 
and d-atoms, but not for heavy metals or the lanthanides. 
Even in those groups where the GGA improves upon the 
LSDA results systematically, there still remains a signif- 
icant difference with respect to the experimental data. 



C. The Spin. 

The current section presents the spin values obtained 
ib-initio for various atoms and ions. For all the systems 
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FIG. 2: (Color online). The first ionization energy as a function of the atomic number Z, for the LSDA (rhombi) and the GGA 
(circles) calculations, compared to the experiment (squares). 
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FIG. 3: (Color online). The relative error in the ionization energy by the LSDA and the GGA calculations comparing to the 
experiment, as a function of the atomic number Z. 



from the s- and p-blocks the calculated spin reproduces 
the experimental value. This result supports the well- 
known empirical Hund's rule 37j that states that the 
relative orientations of the spins of the electrons have to 
be determined to maximize the z-projection of the total 



spin of the atom, S. For atoms from the s-block, S 



,0 



1, §, 1, i in each 



in each row, and for the p-block S - 
row. 

All atoms and ions in the s- and p-blocks, except 
Ba + , are non-interacting pure-state w-representable (NI- 
PSVR; sec 30] for definition), namely, their Kohn-Sham 
electronic configurations are identical to the empirical 
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ones (l €31 ] - 

The calculated and experimental spin values for tran- 
sition metals and their first ions are depicted in Fig. |4] 
From the figure one can see that the number of cases 
for which the spin is predicted incorrectly in the LSDA 
increases with the atomic number. That is, for the first 
row one finds 2 such cases (Ti and V), for the second 
row - 4 cases (Zr, Rh, Y + and Tc + ), for the third row 
- 10 cases. In contrast to the situation in the s- and p- 
blocks, several atoms in the d-block are NI-EVR, with 
fractional occupation of the KS orbitals. These include 
the well-known examples of Fe and Co, as well as the 
less well-known ions of Sc, Ti, Zr and Hf. Details are 
given elsewhere In addition, we find that Ni, Tc, 

Re and Os are NI-PSVR and their spin is predicted cor- 
rectly, but their electronic configuration does not coincide 
with the empirical one. As already mentioned, the KS 
and the empirical electronic configurations do not have 
to coincide. It is pointed out that for these atoms the 
assumption made in previous calculations [H, [3l[ to use 
the empirical electronic configuration is unjustified. 

Furthermore, for Ti, Fe + and Re + the electronic con- 
figuration that is lowest in energy appears not to be 
proper, but only proper in a broad sense (b.s.). As men- 
tioned above, b.s. configurations obey all the required 
restrictions related to a rigorous definition and differ- 
entiability of energy functionals. However, as in most 
cases the energetically favorable electronic configuration 
appears to be proper, these few cases are pointed out. 
As shown below fSec. IIVF]) , for Ti and Rc + this feature 
depends on the specific parametrization chosen for the 
xc- functional (VWN). 

In lanthanides, most atoms and ions are NI-EVR, ex- 
cept Eu, Yb and their ions, in contrast to other blocks. 
In addition, there are two b.s. -proper cases - Ce and Eu + . 
The ion La + has three degenerate KS levels, 4f^,5rf^ and 
6s^ - a unique example in the periodic table [49J . 

Despite the poor correspondence of the ionization en- 
ergies to the experiment in lanthanides (see Sec. IIVB[) . 
for most of them the spin was obtained correctly. From 
Fig. [5] we see that among the atoms only for Gd the cal- 
culation incorrectly predicted S = 3 instead of S = 4. 
For ions, however, there is a systematic error in the spin 
values for Gd + - Tm + - the calculated spin is lower by 
one unit compared to the experiment. 



D. The Effect of Choosing the Occupation 
Numbers ab-initio. 



Throughout the periodic table, the KS electronic con- 
figuration that was obtained ab-initio in the current 
work, using the AEVRA algorithm, differed from the em- 
pirical one for 27 atoms and 26 first ions, mainly in the d- 
and /-blocks of the periodic table. In total, the ionization 
energy, which depends on both the energy of the atom 
and the ion, is affected in 31 cases In these systems 



the assumption made in Ref. 25J that the empirical oc 



cupation numbers can be taken as input to atomic DFT 
calculations is unjustified, and our findings provide the 
desired correction for the LSDA results of Ref. [25| . Since 
for all other atoms the results in the two studies agree 
within the numerical error (see Sec. |V|, the aforemen- 
tioned cases are compared below to examine the effect of 
an accurate treatment of occupation numbers in the KS 
system on the total and ionization energies, and the spin. 

The total energy for all the affected systems was found 
to be lower than in Ref. [25] on average by a sev- 
eral tens of mhartree, significantly above the numeri- 
cal convergence of 2 /ihartree. A maximal difference of 
140 mhartree for Gd + (Z = 64) and a minimal difference 
of 0.5 mhartree for Re + (Z = 75) were obtained. It is 
concluded that in these systems the minimum in energy 
is obtained with the ab-initio occupation numbers, rather 
than with the empirical electronic configuration. 

As for the ionization energies, some of them were ob- 
tained higher, and some of them - lower than those in 
Ref. (25j |. with a maximal difference of 70 mhartree for 
Pt (Z — 78). The relative errors A re i of the ionization 
energy are plotted in Fig. [6] for the relevant atoms. It is 
evident from the figure that for most transition metals 
an accurate treatment of occupation numbers improves 
the results for the ionization energies (note, however Rh, 
Z = 45). In lanthanides, a new problem in the LSDA 
description of atomic systems, which was 'hidden' be- 
fore by an incorrect assumption on the KS occupation 
numbers, is revealed (see also Fig. [71) . While the ioniza- 
tion energies provided in Ref. [251 ] for lanthanides mimic 
the experimental trend with an accuracy of 3%, in the 
current work the calculated ionization energy drops with 
Z, until reaching Eu (Z = 63), for which it experiences a 
sharp rise typical for half-filling a shell. After a small rise 
for Gd, I(Z) starts decreasing again and rises back only 
for Yb (Z = 70), with both 4/ and 6s shells completely 
full. The relative errors in the ionization energies for the 
lanthanides can reach 17%. As mentioned in Sec. IIVE1 
using the GGA does not improve the results for the lan- 
thanides. 

Finally, applying the empirical electronic configuration 
to the KS system always produces the correct spin, by 
definition. However, this is sometimes at the price of 
leaving the KS system in an excited state, which has 
been shown to be forbidden theoretically [30]. Using the 
AEVRA algorithm, which assures that only EVR densi- 
ties are used in the KS system, reveals many cases among 
the transition metals and lanthanides, where the LSDA, 
as well as the GGA, fail to predict the corre ct sp in of the 
system (see Figs. H 03 Sec. [EVE] below and [if). 



E. The Effect of the PBE-GGA Approximation. 

In the current section we compare numerical results 
achieved with the PBE-GGA [11] versus the LSDA. Let 
us begin with the ionization energies. As could already 
be seen from table |TT] and Fig. [3] PBE-GGA provides a 
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FIG. 4: (Color online). The spin S for neutral atoms (left) and first ions (right) throughout the d- block, calculated within the 
LSD A (squares), compared to the experimental data (rhombi). 



systematic, but modest, improvement over the LSDA. A 
major improvement is noted for Fe and Co (Z — 26, 27), 
for which the relative errors drop from 4% and 3% to 
0.7% and 0.2%, respectively, and for Nb-Ru (Z = 41 - 
44), for which the average errors drop from 3% to 0.7%. 



The improvement of GGA for ionization energies is 
illustrated in Fig. [8] For this purpose define R = 
A^ A /A^ DA as the ratio between the GGA and the 
LSDA relative errors. If R < 1, GGA improves over the 
LSDA; if not - LSDA produces a better result. R is de- 
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FIG. 5: (Color online). The spin S for lanthanide neutral atoms (left) and their first ions (right), calculated within the LSDA 
(squares), compared to the experimental data (rhombi). 
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FIG. 6: (Color online). The relative errors A re i in the ion- 
ization energies, with respect to the experiment, for selected 
atoms (see in text), as obtained in the current work and in 
the study by Kotochigova et al. [2o| . 



FIG. 7: (Color online). The ionization energies for the lan- 
thanides (Z — 57 — 70) as obtained in the current work 
(rhombi), comparing to calculations by Kotochigova et al. [25l ] 
and to the experiment [161 ]. 



picted on a logarithmic scale, as a function of Z. It can be 
seen that for lighter atoms with 1 ^ Z ^ 56, in the vast 
majority of cases GGA improves over LSDA. There are, 
however, 2 s-atoms, 7 p-atoms and 7 <i-atoms for which 
a worse result is obtained. For heavier atoms, especially 
the lanthanides, there are many more cases where the 
GGA fails. 

An analysis of how well does the GGA reproduce qual- 
itative trends of the experimental function I(Z), com- 
pared to the LSDA, produces a somewhat disappointing 
conclusion: from Fig.[2]it can be clearly seen that in many 
cases where the LSDA does not succeed to mimic the ex- 
perimental tendency, the GGA fails as well. It happens 
when the calculated results overestimate the value of / for 
Z = 23, 24 and for Z = 28, 29. The same is true for the 
peak at Z = 46 and a wrong tendency for Z = 51,52. 
Most important, for the lanthanides the GGA has the 



same shortcomings as the LSDA. For atoms with Z > 70, 
for which a relativistic treatment is required, the GGA 
provides results qualitatively similar to the LSDA, with 
large discrepancies comparing to the experiment. 

Comparing the spins and the electronic configurations 
produced by the GGA, with respect to the LSDA, we 
find little difference. In the d- and /-blocks, only in two 
cases a different spin value is obtained within the GGA: 
for the ion Re + GGA correctly reproduces S = 3 and 
for the ion Gd + the GGA predicts S = |, while the ex- 
perimental result is | . For other cases wrongly obtained 
in the LSDA, no improvement in the prediction of the 
spin occurred. Almost all atoms and ions which are NI- 
EVR in the LSDA, are NI-EVR also in the GGA (except 
Gd+ and Dy+, which are NI-PSVR in GGA). In addi- 
tion, the ion Re + and the Ti atom, which were obtained 
b.s. -proper in LSDA, are now fully proper. The ions Tc + 
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FIG. 8: (Color online). The ratio R (denned in text) vs. the 
atomic number Z, on a logarithmic scale. Rhombi correspond 
to atoms with R < 1, squares - to atoms with R ^ 1. 
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FIG. 9: (Color online). The differences of the ionization 
energies obtained with LSDA(PZ) (rhombi), LSDA(PW92) 
(squares), GGA(VWN) (x's), GGA(PZ) (hollow rhombi) 
and GGA(PW92) (hollow squares), with respect to the 
LSDA(VWN) results. 



and Gd + become b.s. -proper in the GGA being proper for 
LSDA. The ions Fe + , Eu + and the Ce atom remain b.s.- 
proper. The fractional occupation numbers for NI-EVR 
cases are somewhat different, of course. This observation 
supports the statement [13] that emergence of NI-EVR 
densities is not an artifact of the LSDA approximation, 
as was claimed in the past [20j . Detailed tables with the 
spin values and electronic configurations are available on 
EPAPS go[. 



F. Other Parametrizations to the LSDA 
Correlation Energy. 

Several analytical expressions for the correlation en- 
ergy functional have been derived and fitted to the 
Monte-Carlo calculation results for the homogeneous 
electron gas obtained by Ceperley and Alder [Pol ]. 
Three expressions are popular in the literature: the 
parametrization by Vosko, Wilk and Nusair (VWN) [Toj . 
by Perdew and Zunger (PZ) [ll[ and by Perdew and 
Wang (PW92) [l2[ ■ The results presented so far all used 
the VWN parametrization. In the current section numer- 
ical results obtained with all three parametrizations are 
compared, within both the LSDA and the PBE-GGA, for 
atoms with 1 ^ Z ^ 56 and 71 ^ Z ^ 86, excluding the 
lanthanides. Since all three parametrizations are based 
on the same Monte-Carlo calculation [lEj], their results 
are expected to be very close. 

Figure [9] presents the differences of the ionization en- 
ergies with respect to the LSDA(VWN) results: AI = 
jx _jLsda(vwn) _ Indee d, one can see that all the LSDA 
results differ very little from each other. The PW92 is 
much closer to the VWN, with a maximal discrepancy of 
0.4 mhartree for He. The PZ has larger discrepancies, 
but most results lie within a scatter of 1 mhartree. A 
very similar picture is obtained within the GGA results. 



The question of whether the discrepancy in I(Z) ob- 
tained with different parametrizations overshadows the 
difference between the LSDA and the GGA is answered 
positively. In Fig. [9] it is clearly seen that the GGA re- 
sults remain distinct from the LSDA results for almost 
all atoms, except maybe Z = 7, 11 — 15, 33. 

As for the spin, the only affected system was Fe + , 
for which LSDA(PZ) produced S = §, while other 
parametrizations gave S = i. In addition, Ti, 
which was obtained as b.s. -proper for LSDA(VWN), is 
proper for all other parametrizations, as well as for 
the GGA parametrizations. Re + was obtained b.s.- 
proper for LSDA(VWN) and LSDA(PW92), but proper 
for LSDA(PZ) and all the GGA parametrizations. The 
Os atom, on the contrary, appeared b.s. -proper for 
LSDA(PZ) and proper for all the rest. No connection 
between these isolated cases was found. 



V. DISCUSSION AND SUMMARY. 

In the present contribution, self-consistent ab-initio 
non-relativistic DFT calculations for atoms and ions with 
1 ^ Z ^ 86 have been reported. The total and ionization 
energies, the spin and the Kohn-Sham electronic configu- 
ration were obtained within both LSDA and PBE-GGA. 
The calculations were performed in three parametriza- 
tions for the correlation energy [nj{l2j ■ A numerical con- 
vergence of 2 /mhartree in the total energy was achieved. 

Comparing the LSDA results presented here with pre- 
vious calculations reported in Ref. [25| . we can distin- 
guish between two cases. In cases where the electronic 
configurations are the same, we find that the discrepancy 
between the total energies does not exceed 2.5 /mhartree, 
and between the ionization energies - 3 /ihartree. These 
findings are in agreement with the numerical accuracies 
reported in both studies. However, in cases where the 
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configurations are not the same, the discrepancy in the 
total energies, and therefore in the ionization energies, 
can be much larger, reaching 70 mhartree. This differ- 
ence is related to the restrictive nature of the energy min- 
imization undertaken in Ref. [25| . which is overcome in 
the present work by employing the recently proposed [30] 
algorithm AEVRA. The analysis conducted in Sec. IIVDI 
showed that the assumption that the KS electronic con- 
figuration is identical to the empirical one [25[ is unjus- 
tified for many systems from the d- and /-blocks of the 
periodic table. Once this assumption is lifted, calculated 
ionization energies for several transition metals fit the 
experimental values better. 

The GGA results of this work can be compared to 
Ref. [HI, which used the PBE approximation [l4| to ob- 
tain results for atoms in the s- and p-blocks. For all 
atoms the differences in total and ionization energies are 
smaller than the numerical convergence error of Ref. [26[ . 

The shortcomings of the currently used approxima- 
tions to the xc-energy functional are markers against 
which future improved xc-approximations can be mea- 
sured. In this context, we discuss the numerical results 
presented in the current contribution. 

In Sec. lIVBl we saw that within both the LSDA and the 
GGA ionization energies are reproduced for many atoms 
within a few percent. However, for atoms with Z > 70 
the average accuracy drops to about 10%. In addition, 
the calculated ionization energies in the lanthanides show 
a qualitatively different trend from the experimental re- 
sults, with a maximal discrepancy of 13%. Poor results 
are obtained also for several transition metals, e.g. V, Cr, 
Ni, Cu, Rh and Pd. 

In Sec. lIV Cl we showed that the spin of atoms and ions 
was correctly reproduced throughout the s- and p-blocks. 
Several mismatches with respect to the experimental spin 
values occurred in the (i-block, especially for the heavier 
atoms. For the lanthanide ions with Z = 65 — 69, the 
calculated spin was systematically lower than the exper- 
imental one by one unit. 

In Sec. IIVEI an analysis of how PBE-GGA affects the 
results for atomic systems was presented. It was found 
that although the PBE-GGA slightly improves ionization 
energies of atoms, it fails to repair the main shortcomings 
of the LSDA. In addition, the spins of atoms and ions 
remained essentially unaffected. 

In Sec. IIVFI results obtained using each of the three 
parametrizations of the correlation energy - VWN [l(|, 
PZ [TO] and PW92 [H[, both within the LSDA and the 
GGA, were compared. The discrepancy in the ionization 
energies obtained with different parametrizations was 
found to be small, of the order of 1 mhartree. This dis- 
crepancy provides an error estimate of the parametriza- 
tion process, and as a result, a lower bound for the 
error of the whole calculational scheme. The spin of 
atoms and ions was almost insensitive to the choice of 
the parametrization. 

Several reasons can be given for the observed discrep- 
ancies in the ionization energies and spins of atoms and 



ions. First, relativistic effects, which could not be repro- 
duced in the current treatment, are relevant for heavy 
atoms. For example, the dip in the ionization energy that 
is typical for the 4 th atom in each of the first four rows of 
the p-block (see Fig. [2]), occurs for the 3 rd , rather than 
the 4 th , atom (Bi, Z = 83) in the fifth row. This can be 
explained considering the j-j coupling in heavy atoms, 
which cannot be reproduced in non-relativistic calcula- 
tions. Furthermore, in the third row of the d-block the 
calculated I(Z) presents sharp increases for Z = 78, 80 
and a drop for Z = 79, which can be related to filling 
the (i-orbital at Z — 78 and then the s-orbital at Z = 80. 
Relevant for the second row of transition metals, this be- 
havior is not observed experimentally in the third row, 
suggesting it has to be treated relativistically, similarly 
to the fifth row in the p-block. It can be also observed 
that the calculation provides the same spin values for the 
corresponding atoms and ions from the second and third 
rows of the d-block. In the second row the correspon- 
dence with the experiment is high, which is not the case 
in the third row. It is probable that the relativistic ef- 
fects in the third-row atoms are strong enough to dictate 
a different spin value. 

Second, the spherical approximation to the density em- 
ployed in the current calculations can affect the results 
in the p-,d- and /-blocks. Without such an approxima- 
tion, eigenenergies £ n ima arid eigenfunctions ipnima with 
different axial numbers m — would become dis- 

tinct. The manifolds of these states may then overlap, 
which will affect the occupation of the levels and also the 
density n(r). In particular, systems that in the spherical 
approximation are NI-EVR, may become NI-PSVR. In 
addition, the angular momentum L of the atom can be 
obtained and compared to the experiment. Furthermore, 
in cases where the gap between the occupied and the un- 
occupied orbitals with opposite spins is small, the spread 
between orbitals with different m-numbers can cause a 
change in the spin of the atom. This can possibly cor- 
rect the prediction of the spin values for some d- and 
/-systems. However, early studies [42|, |43| suggests that 
a non-spherical treatment with local xc-functionals will 
only slightly affect the physical results. Furthermore, it is 
possible that calculations beyond the spherical approxi- 
mation would require much better approximations to the 
exchange-correlation functionals - the so-called orbital- 
dependent functionals [44j . but at present they are com- 
putationally more expensive [45J and beyond the scope of 
the present report. Moreover, while the exact exchange 
potential only is available, an accompanying correlation 
potential is still under development. 

Third, the local or semi-local nature of the xc-energy 
employed may lead to additional inaccuracies. One of the 
manifestations of (semi-)locality is the spurious interac- 
tion between an electron and itself, which can influence 
the results all over the periodic table. 

In particular, the incorrectly predicted spin for the 
lanthanide ions, as well as the fact that NI-EVR solu- 
tions occur in lanthanides frequently, raise the suspicion 
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that the 6s' level is obtained with a too high eigenen- 
ergy with respect to the 4/ levels, probably due to the 
self-interaction of electrons. 

To investigate this hypothesis, numerical checks have 
been performed on selected lanthanides. To lower the 
eigenenergy e 6s t artificially, the nuclear charge was raised 
from Z to Z + S. The value of S was chosen to be the 
lowest possible to separate the degenerate Ap and 6s^ 
energy levels, both in the atoms and in the correspond- 
ing ions, and thus produce NI-PSVR results. In Pm, 
Nd and Sm (Z = 59,60,62), 5 was 0.180,0.301,0.094, 
respectively. The ionization energies obtained are / = 
0.204, 0.228 and 0.187hartree, correspondingly, which are 
considerably higher than the original results, and thus 
closer to the experimental values. However, the spins of 
these atoms are now obtained incorrectly, all larger by 
one unit with respect to the experiment. The reason is 
that the formerly occupied 6s^ level becomes unoccupied 
due to the increase in Z. Furthermore, a similar check 
for Tm (Z = 69, 8 = 0.35) raised the ionization energy 
to 0.255 hartree, but did not affect the spin of its ion, 
which was obtained incorrectly. These tests have only 
a qualitative character, showing that proper treatment 
of the self-interaction of electrons might improve the re- 
sults for I(Z) in lanthanides, but improvement in S is 
not assured. 

In addition, the question of how stable is a calculated 
spin value for a given atom/ion with respect to small 
variations of the xc-approximation can be addressed us- 
ing the current database. To do so, we define a criterion 
of stability 

D(S) = mhi(£(S') - E(S)) (17) 

as the difference between the ground-state energy E(S) 
with the nominal spin S, and the ground-state energy 
of the same system given it has the spin S', which is 



the closest to E(S). D is an estimate of the magnitude 
of variation in the total energy differences required to 
change the spin of the system. 

Examining the d- and /-blocks, it was found that for 
open-shell atoms D can reach values of 0.08 — 0.10 hartree 
for e.g. Cr + , Tm, Mo + , Ni + , as well as values of 
0.001 hartree and lower for Fe + , Gd+ and Re + . No 
correlation was found, however, between systems with 
a low D and an incorrect prediction of the spin or a 
high relative error in the ionization energy. For some 
atoms with an incorrect spin the result is very stable, e.g. 
D(S) = 0.07 hartree for Pt and 0.03 hartree for V. There- 
fore, the mismatches in spin in these atoms are probably 
not subject to nuances of the GGA or the parametriza- 
tion of the correlation energy, but must be corrected with 
a rather different approximation to the xc-energy. 

To conclude, we have calculated the energies and spins 
of atoms and ions within the local spin-density approx- 
imation and the semi-local GGA. We find that moving 
from a local to a semi-local approximation to the xc- 
energy does not improve significantly the description of 
atomic systems in DFT. In particular, the spins of some 
transition and lanthanide elements are obtained incor- 
rectly. Further improvement in the description of atomic 
systems in DFT should focus on non-spherical calcu- 
lations and on advanced non-local exchange-correlation 
functionals. 
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